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Abstract 

Landau theory and symmetry considerations lead us to propose an explanation for several seem- 
ingly paradoxical behaviors of charge ordering (CO) and spin ordering (SO) in the mixed valence 
compound LuFe204. Both SO and CO are highly frustrated. We analyze a lattice gas model of CO 
within mean field theory and determine the magnitude of several of the phenomenological interac- 
tions. We show that the assumption of a continuous phase transitions at which CO or SO develops 

> 

implies that both CO and SO are incommensurate. To explain how ferroelectric fluctuations in the 

o ' 

charge disordered phase can be consistent with an aniiferroelectric ordered phase, we invoke an 

(N 

■ electron-phonon interaction in which a low energy (20 meV) zone-center transverse phonon plays a 

o, 

key role. The energies of all the zone center phonons are calculated from first principles. We give a 
Landau analysis which explains SO and we discuss a model of interactions which stabilizes the SO 
^ . state, if it assumed commensurate. However, we suggest a high resolution experimental determi- 

nation to see whether this phase is really commensurate, as believed up to now. The applicability 
of representation analysis is discussed. A tentative explanation for the sensitivity of the CO state 
to an applied magnetic field in field-cooled experiments is given. 

PACS numbers: 75.25.+z,75.10.Jm,75.40.Gb 
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I. INTRODUCTION 



The phenomenon of charge ordering (CO) has been studied ever since the observation 
of the Verwey transition jlj in Fe^C^ (in which the average valence of the Fe ions is 8/3). 
An oft-cited paper by Anderson [2j proposed a simple and appealing model to explain the 
ferroelectric charge ordering. However, the most recent high-precision neutron scattering 
resultsjsj show that this model does not correctly explain the CO in Fea04 and the expla- 
nation of the actual nature of its CO remains an open question.^] Similarly i?Fe204, where 
R is a trivalent rare-earth and in which the average valence of the Fe +2 and Fe 3+ ions is 5/2 
presents an even more challenging problem to our understanding of CO. In this paper we 
will consider LuFe 2 4 (LFO) whose CO and magnetic structure has been widely studied in 
recent years. (5- 13]. 

As an introduction we review the most salient experimental results relevant to LFO. In 
Fig. [1] we show the trigonal lattice structure i?3m[14|Jl5| of LFO. The generators of the 
point group are a) inversion about the center of the unit cell, b) the x-y mirror plane, and 
c) the three-fold axis. Note that the Fe ions form triangular lattice layers (TLL's) arranged 
in bilayers. The bilayers are separated by a TLL of Lu ions. The stacking of the Fe TLL's 
is in the same order as for an fee crystal. The rhombohedral unit cell spans three bilayers 
and there are two Fe ions per rhombohedral primitive unit cell, as shown in Fig. 1. At 
temperatures above 500K, the valence electrons can thermally hop so that effectively all Fe 



ity, the dielectric 

3- 



sites appear to have charge 2.5e. 16] Consistent with this electronic mobi 
constant at zero frequency (shown in Fig. [2]) is very large for T > 200 K 

As the temperature is reduced from 500K, CO correlations develop at wave vectors which 
nearly coincide with "root 3" (R3) ordering (see Figs. 5 and 6, below) within each TLL 
and eventually at the charge ordering temperature Tqq 320K three dimensional long 
range CO develops via a continuous transition. |l3|] Both the fluctuations and the long range 
order occur at incommensurate values of the wave vector, [sjl^ In the paraelectric phase 
(T > Tco) the dominant fluctuations are consistent with no enlargement of the unit cell 
in the c direction. 12J We call such fluctuations "ferro-incommensurate" (FI) fluctuations 
to emphasize that their wave vector has incommensurate in-plane components. (The in- 
commensurate in-plane components are very close to the values of the X point of a two 
dimensional triangular lattice gas with repulsive interactions, as we discuss below.) Surpris- 
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FIG. 1: (Color online) Left: Fe ions in LFO. The brackets indicate a bilayer consisting of two Fe 
RLL's. The presence of a TLL of Lu ions between adjacent Fe bilayers is indicated. The oxygen 
ions (not shown) are almost uniformly distributed over the structure. The hexagonal (conventional) 
unit cell contains three bilayers configured so that the TLL's are stacked in the order ABCABC 
(A=red, B=black, C=blue) with two sites in the rhombohedral unit cell labeled "1" and "2" which 
are related by a center of inversion symmetry. The screened Coulomb interactions U n used below 
in our calculations are indicated, d = 2.52A and d! = 5.81A 



ingly, the CO that occurs for T < Tqq involves a doubling of the unit cell along c. |12lll7 ] 
We call this ordering "antiferro-incommensurate" (AFI). These structures are shown in Fig. 
[3] where, for simplicity, the small incommensurability of the wave vectors is neglected. A 
main objective of the present paper is to explain why the CO phase is AFI and does not 
reflect the dominant FI fluctuations of the paraelectric (P) phase. 

We will analyze this unusual CO within the lattice gas model used by Yamada et a/[6|] 
which we refer to as Y. The most striking result found by Y was that even with in-plane 
coupling U\ and interplane couplings U% and U 3 , long-range order is not possible because the 
maximum of the wave vector dependent susceptibility occurs over an entire 'degeneration 
line' in wave vector space. This result had been known for similar spin and lattice gas 
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FIG. 2: The dielectric constant at zero frequency from Subramanian et a/ 9|] 
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FIG. 3: (Color online) From Ref. 



12J. Fully CO structures when the incommensurability of the wave 



vector is neglected. The dipole moments of each bilayer (in the absence of incommensurability) as 



calculated in Ref. 
R3 structure. 



12j are shown. Left: The 'antiferroelectric' R3 structure. Right the "ferroelectric' 
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191 ] and later of 



models on a rhombohedral lattice from the work of Rastelli and Tassi 
Reimers and Dahn[20|. As Y found, it was necessary to include an interaction (U4 in Fig. 
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1) between next-nearest neighboring TLL's in order to remove this degeneracy. We will 
analyze this situation in detail and show that there are two crucial parameters which govern 
this phenomenon. The first parameter is the interaction £7 2 in Fig. [T] which scales the radius 
of the cylinder on which the degeneration line is wrapped. The second parameter is the 
interaction U4 in Fig. [TJ which removes the continuous degeneracy and leads to an energy 
variation of amplitude AE as one traverses the degeneration line. To elucidate our mean- 
field results we will briefly review the results for a single TLL with repulsive interactions. 
Repulsive interactions are clearly relevant for the charge-charge interactions. For magnetic 
interactions the presence of many exchange paths through intervening oxygen ions suggests 
that the magnetic interactions should also be repulsive (i. e. antiferromagnetic). The 
most general result of our analysis is that the wave vector of the stacked TLL's of LFO 
is unstable at X (the wave vector which characterizes the R3 structure) if a continuous 
transition is assumed, in which case the ordered phase must perforce be incommensurate. 

This same logic applies to the magnetic phase transition into the spin ordering (SO) phase 
which appears at T = 240K. We will discuss the ramifications of the fact that the magnetic 
transition appears to be a continuous one to a commensurately ordered spin state. We 
will give a Landau analysis of the symmetry of the SO phase and will discuss microscopic 
interactions which can explain this ordering. Finally, we will discuss briefly a possible 



131 ] which show that such a protocol has seems 



explanation of recent field-cooled experiments 
to significantly destabilize the AFI CO state. 

Below 170K the system undergoes another transition in which the magnetic order pa- 



rameter sharply decreases. [Hi] The details of this state are not settled at present |l3| and we 



will not consider it further. Thus the magnetoelectric phase diagram we are considering is 



that shown in Fig. HI A brief summary of this work appeared some time ago. 21] 



II. MEAN FIELD TREATMENT OF CO 



Calculation 



We start with a Landau analysis of CO using the lattice gas model of Y in which one 
introduces a variable a n (R), where n (n = 1, 2) labels the nth site of the rhombohedral unit 
cell at R. which assumes the value +1 (—1) if the site is occupied by an Fe 3+ (Fe +2 ) ion. 
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FIG. 4: The magnetoelectric phase diagram of LFO based on the results discussed in the text. 

Then x n (R) = (cr„(R)), where ( ) is a thermal average. As shown in Fig. 1, we include an 
interaction U\ between nearest neighbors within the same TLL, an interaction TJi between 
nearest neighbors in adjacent TLL's within the same bilayer, an interaction {7 3 between 
nearest neighbors in adjacent bilayers, and an interaction {7 4 between nearest neighbors in 
second-neighboring TLL's. As argued by Y, in view of the large dielectric constant (See 



Fig. [2]) we prefer to use a model in which the interactions : 
rather than one based on a long-ranged Coulomb interaction 
in terms of the Fourier transformed variables 



al 



off rapidly with separation 



221 ]. The free energy is written 



£n(q) = ^2x n (R)exp(iq-Ii) . (1) 

R 



Then 



F = \Y. E ^ m (qK(q)*x m (q) + 0(x 4 ) , (2) 

q n,m=l 

where F nm (q) = F mn (q)*. The free energy must be invariant under spatial inversion X, since 
X is a symmetry of this lattice. Under spatial inversion X, site #1 goes into site #2, so that 
XS'i(q) = 52 (q)*, from which we conclude that Fu = F 2 2- A continuous CO transition is 
signaled by the appearance of a zero eigenvalue of the quadratic form. This instability will 
first occur at a wave vector whose value we wish to determine. 

Before analyzing the model in detail we review some results for the simpler problem of a 
lattice gas with repulsive interactions on a single TLL. In left panel of Fig. |5]we show the 
first Brillouin zone for the TLL with the X points which are the wave vectors of the ordered 
phase of this system. Note that the X point (which gives rise to the CO or SO structure, 
shown in the right panel of Fig. [5]), is an isolated point having higher symmetry than that 



X = (1/3, 1/3) RX = X 



FIG. 5: (Color online) Left: The reciprocal lattice of the TLL. The X points all equivalent to one 
another under a three- fold rotation R, as discussed in the text. Right: the R3 structure associated 
with the wave vector of the X point. The amplitudes Z of the three sites within the R3 unit cell 
are given by Za = Z cos((j)), Zb = Z cos((f> + 2n/3), and Zc = Z cos(4> + 4tt/3). The choice <p = 
yields {Za, Zb, Zc) oc (2, —1, —1) and the choice <f> = ir/6 yields (Za, Zb, Zc) oc (1, 0, —1). 



of all surrounding points. To see this, note that uniquely this point is invariant under a 
three-fold rotation about the center of zone, because under the three-fold rotation such a 
point is taken into another X point which is equal (modulo a reciprocal lattice vector) to to 
the original point. As a consequence of the special symmetry of the X point, its wave vector, 
determined by the instability in the quadratic term of the free energy, is stable with respect 



to the addition of small further neighbor interactions. |23| Such structures have been observed 
many times over the last half century. 24| In contrast, consider the analogous X points in the 
rhombohedral reciprocal lattice shown in Fig. [6j (We label these points as X regardless of 
the value of q z and refer to them as an "X line.") Here, an X point is not invariant under a 
threefold rotation because the points before and after such a rotation are not equal modulo 
a reciprocal lattice vector. (The point is that the rhombohedral reciprocal lattice vector 
does not connect points before and after a three-fold rotation because the reciprocal lattice 
vector needed to relate the components in the plane does not have q z = 0.) So if a transition 
at this wave vector is continuous, the wave vector must perforce be incommensurate. Thus, 
without any calculation, we have shown that Landau theory explains why the CO wave 



■ q z = • q z = 271 /c q z = - 271 /c 



X = (1/3, 1/3) R X = X' 

G X' = X + (2 7C/c)k 

FIG. 6: (Color online) The reciprocal lattice of the rhombohedral lattice. 



vector is incommensurate. (How this conclusion relates to representation theory is discussed 
in Appendix A.) Accordingly, the CO phase should not display a spontaneous polarization, 
P. The nonzero value of P may be an artifact of the small electric field applied during 
the experiment (as argued in Ref. Il2l ). Also, it is possible that the pyroelectric current 
(from which the value of P is deduced) might be confused with currents which develop 
in the conductive sample. 25[ At a recent conference it was reported that P = at low 
temperature. We will discuss below that a similar analysis of the magnetic phase transition 
at T = 240K indicates that if that transition is continuous the ordered phase ought to be 
incommensurate. 

We now return to the explicit calculation of the incommensurate CO wave vector. As 
mentioned in the introduction, from previous work it is known that the minimal 

model that gives stable three-dimensional long-range order requires the interactions shown 
in Fig. 1. Other interactions (such as a second- neighbor in-plane interaction or a second 
neighbor interaction between adjacent TLL's) only lead to perturbative corrections. There- 
fore we simplify the analysis by only considering the minimal model. To determine the 
wave vector of CO within this model, we analyze the appearance of a zero eigenvalue of the 
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quadratic form of Eq. (j2]) . In mean field theory one writes 

F nm (q) = c'kT5 nm + ^ U(0, m; R, n) exp(zq • R) , (3) 

R 

where d is a constant of order unity and U(0, m; R, n) is the interaction between sites m 
in the rhombohedral unit cell at the origin and n in the rhombohedral unit cell at R. We 
set d — 1, and fee = 1, and henceforth, unless stated otherwise, all energies will be in 
temperature units. In Cartesian coordinates 

F\\ = T + Ui[2 cos(aq x ) + 4 cos(aq x /2) cos(y/3aq y /2)] 

+U 4 [e lc ^ 3 A(q x , q y ) + e- ic ^ 3 A(q x , q y )*} 
F 21 = U 2 e- 2lc ^ 3 A(q x , q y ) + U 3 e^ 3 A{q X) q y )* , (4) 



where 



A(q x ,q y ) = 2e iq y a ^ /6 cos(aq x /2) + e^ 1111 ^ 3 . (5) 



To organize the calculation we will consider U n /U\ for n > 1 to be of order the expansion 
parameter A and we will work consistently to the lowest sensible order in A, keeping in mind 
that for (74 = we have a line of infinite degeneration^ 
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20| | . In any event these works 



indicate that for U\ > the instability in the quadratic form first appears near the X point, 
for some discrete values of q z , for ?7 4 ^ and for all q z for U± = 0. Accordingly, we set 
aq x = Ait / 3 + p x and aq y = p y and determine p to leading order in A for arbitrary q z . 
For this purpose we evaluate the matrix F nm . We find that up to order £7iA 2 , 

F n = T + Ui[-3 + (3/4)p 2 ] - ^ 4 v / 3pcos(0 + ck z /3) , (6) 

where pcoscj) = p x and psin0 = p y , with p > 0, and p 2 = p 2 x + p 2 . As we shall see below, 
p = 0(A) and to clarify the situation it is only necessary to calculate the eigenvalues of the 
quadratic form F to order X 2 Ui. 

A = -( v / 3/2)pe^ . (7) 

Thus the critical eigenvalue, p(q) which first approaches zero as the transition is ap- 
proached is given by F u — |F 12 | or, up to order Ui\ 2 , 

p(q) = r + t/ 1 [-3 + (3/4)p 2 ]-t/ 4 v / 3pcos(0 + cg,/3) 

-^TT \ U l + U l + 2U ^ cos ( 2< ^ - c ^/ 3 )l V2 • ( 8 ) 



Note that for a single TLL, for which U 2 = U3 = U4 = the X point with p = is stable 
and that, in view of the term linear in p, the X point becomes unstable in the presence of 
interlayer interactions, One might have expected to have three dimensional long-range 
order when Ui, U 2 , and U 3 are all nonzero because then each TLL's interacts with ones 
above and below it. However, the special symmetry of the rhombohedral lattice prevents 
when only these interactions are present, so we are forced to include a 



ordering [6 



nonzero value of U4. 

We first minimize /i(q) with respect to p, which is determined by 

= % = ~^ Ua cos(0 + cqz/3) + \ Ulp ~ T'P* + u " + 2U2U3 cos(20 ~ cqz/3)] • (9) 

Thus the value of p which minimizes p and which we denote p* is given by 



p* = (2v / 3/3)(t/ 4 M)cos(0 + cg 2 /3) 



^/3 

PI + Ul + 2U 2 U 3 cos(20 - cq z /3)} 



3fA 



(10) 



As we mentioned, p* becomes nonzero at order A. Since dp /dp = at the extremum, 
corrections to p* at the next order in A do not affect the result we find for the critical 
eigenvalue p(q z , (ft, p*(q z , 4>))\ 



p{q z ,<f>) = T-3U 1 -V3p*U4Cos(<j) + cq z /3) + (3/4)U 1 (p*) 2 
-{V3p*/2)[Ul + Ul + 2U 2 U 3 cos(20 - cq z /3)} 1 ' 2 
= T — 3Ui — (Uf/Ui) cos 2 (0 + ) - (U 4 X/Ut) cos(0+) - X 2 /^ 



where 



+ ck z /3 and 



X 



u 2 + ul 



2U 2 U 3 cos(20 - cq z /3) 



1/2 



(12) 



When U4 = 0, the critical eigenvalue is a function of the variable 20— cq z /3 and the eigenvalue 
is minimal for 2(fi — cq z /3 = or 2n if U 2 U 3 > and for 2<f> — cq z /3 = n or — n if U 2 U 3 < 0, 
consistent with the results cited for the line of degeneration. 

The values of and cq z which complete the determination of the critical wave vector 
when U4 7^ are selected as those which minimize p. If we define 



n = UJ\U 



41 ) 



m 



u 2 u 3 /\u 2 u 3 



(13) 
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TABLE I: The critical wave vector as a function of the f/'s is given by q = (pcos </>, p sin phi, q z ), 
where p is given by Eq. ()10p and eft and q z are given below, where p is an integer. 



U 2 Ua 
\U 2 U 3 \ 


C/4 




cq z 


+ 


+ 


~~3~ 


—2pir 


+ 




(2p+l)7T 

3 


2vr(-p+ 1) 




+ 


(2p+l> 

3 

2tt(p+1) 
3 


2tt(-p-1/2) 
2ir(-p + l/2) 



then we see that /i is minimized by setting 
cos(0 + cq z /3) = n , 

so that 

+ cg z /3 = (1 -n)7r/2 + 2£;7r , 



cos(20 — cq z /3) = m , 



(14) 



20 - = (1 - m)vr/2 + 2/vr , (15) 



where the integers k and I are free parameters. We therefore get the results of Table [J 
shown in Figs. [TlfTUl 



B. Comparison to Experiment 

Now we compare these results with experiment and if we obtain agreement we should 
be able to identify some of the parameters. Look at the data shown in Fig. [TTJ Note that 
the diffraction pattern at I = 15 (or I = 18) can be compared with that for I = in Figs. 
[7] to [KB Note that the data indicate that the diffraction at I — 15 occurs for <ft = n [since 
it is closer to (0,0,/) than is the commensurate location]. So this is the same as shown in 
Fig. [8] for U2U3 > and U4 < and can also be confirmed by comparison with that case 
(for p — 1) in Table I. Furthermore, as one moves in the direction of positive I = q z , one 
sees that d<j)/dl = —2tt/3 in both the experiment and in Fig. [HJ (Although this line shows 
a definite sign of helicity, the system as a whole is not chiral. Of the six lines equivalent to 
the one shown in Fig. [TlJ three of them have one sign of helicity and three the other sign 
of helicity.) We can also check that our eigenvector of the matrix F of Eq. (j2J) agrees with 



that used in Ref. 
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They use [1,-1] for FI (I = 0) diffraction (see the first paragraph of 



11 




FIG. 7: (Color online) Scattering for positive U2U3 and positive U4. The x and y axes of the 
scattering vector q are indicated. Its z-component is q z = 2irz/c. The X point, which in hexagonal 
notation is (1/3, 1/3,0), is indicated. Here we show the diffuse scattering in the charge disordered 
phase (where the colors have no significance.) If the system were to order, one would have Bragg 
scattering from each of the three domains. The scattering from a single domain is indicated by a 
single color (red, blue, or gold). Note that the diffuse scattering exhibits all the symmetries of the 
crystal, whereas a single domain has lower symmetry, since it does not have three fold symmetry. 
However, if all domains are equally populated, then the three fold symmetry is restored. This is 
the result for the crystal shown in Fig. 1. Note that the mirror which takes y into — y is not a 
symmetry of the crystal and it leads to a slightly different diffraction pattern. 
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FIG. 8: (Color online) As Fig. [7J but for positive U2U3 and negative C/4. 



p3). For us to obtain that result _F 12 must be real positive (to give a minimal eigenvalue). 
This implies that = 71 which agrees with Table I. It is more problematic to connect the 
AFI diffraction to our analysis because the AFI phase can not be explained by the present 
theory, although the AFI diffraction is that shown in Fig. [TU1 

Now we fix the parameters to fit the existing data. In this connection we will assume that 
U\ ^> U2 ^> U3 ^> U4, as will be justified a posteriori. Then the mean-field value of the CO 
transition temperature, which we denote T co ,mf is determined by setting /i(q) = 0. Since 
we believe that the heavily screened interactions decay rapidly with distance, this condition 



13 



FIG. 9: (Color online) As Fig. [7J but for negative U2U3 and positive U4,. 



gives 

Tco,mf = 3C/i . (16) 

If we were to identify this with the observed CO transition at Tqo = 320K, then we would 
conclude that U\ = 110 K. But since the coupling between bilayers is very weak, the two 
dimensional fluctuations will cause the observed value of Tqq to be very much less than 



7co,mf- Accordingly we adopt the estimate [27j 



U x = 500K . (17) 
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z = 



z = - 1 



FIG. 10: (Color online) As Fig. but for negative U2U3 and negative U4. 



Next, we have to decide whether to take U 2 U 3 to be positive, as suggested by the fact 
that the P phase diffraction dominantly occurs at integer values of cq z /(27r) 12j or to be 
negative, as s ugges ted by the fact that the CO phase diffraction occurs at half-integer values 
of cq z / (27r).|s 812] Since the corrections to our mean field theory are smallest in the P phase, 
we use the P phase data to fix U2U3 > and hope to explain the CO data by some correction 
to this theory. (Of course, in addition, if we took U2U3 < 0, we would have to explain why 
screening causes U 2 or U 3 to be negative.) Equation (flUI) gives 



U 2 = VZpUx = O.O6E/1 = 30K . 



(18) 



15 



11(120) 



11(110) 



FIG. 11: (Color online) From Ref. Il2l . Diffraction maxima in the charge disordered phase (T 
360K) near an X line. Hexagonal coordinates are used. Different colors correspond to wave vectors 
of different domains when CO occurs. 



In view of the large ratio U\jU% = 15, it seems reasonable to get U 2 under the assumption 
that U3/U2 <C 1. We do not have an unambiguous way to determine U3 and U4. However, 
considering that U 2 /Ui is about 15, we guess that U^jU^ = 15, which would indicate that 
IC/4I = 2K. Previously we determined that to fit the diffraction at T = 360K we needed to 
assume that U4 was negative, so we set 

U 4 = -2K . (19) 



C. FI versus AFI Transition Temperatures 

Now we want to estimate the difference between the mean-field values of the transition 
temperatures for FI and AFI CO which we denote respectively as Tmf,f and Tmf.af- (By 
our choice of parameters AT = T MF F — T MFj af is positive, which does not agree with the 
experimental value. We have that 

AT = fi(cq z /n) - fi(cq z = 0) , (20) 

when these /x's have each been minimized with respect to and p. Accordingly we need 

0M6g«) 2C/J . U 4 X 

= al = TT~ cos W , + sm W>+ + -J7- sin 0+ 

0(p U\ Ui 
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1 ll 3 [X + 2U 4 cos(0+)] sin(20 - cq z /3) (21) 



2XUt 

For cq z = 0, this minimization gives <ft = 7T, so that 

/i(cg 2 = 0) = T - 3C/! - ZTf/C/i + tf 4 (tf 2 + [7 3 )M - + U 3 ) 2 /{4U X ) . (22) 

Now we analyze the extremum of \x for cq z = ir. All terms in Eq. ( 12 ip are of order X 2 Ui. 
So we assume that U4/U2 and U4/U3 are small and work to first order in those quantities. 
We then find that Eq. (j2ip yields that the extremum occurs for <ft = <p*, where 



7T U4 _ 7T 

6 ~ Th = 6 



(23) 



Then 



M«Z. = vr) = T - 3^ - + 0(f/ 4 2 ) (24) 



and 



Note that the implied negative sign for U4 is crucial to explain the dominance of FI fluctu- 
ations for T > Tqo- Using our admittedly arbitrary estimate of U4 we have 

£^ = _(-2)(30) 
Ut 500 V ; 

In view of the effect of large two dimensional fluctuations we estimate that more realistically 
this model would give 

AT = 0.04ft: . (27) 

Of course, experiment 12| tells us that AT < (i. e. the first criticality we encounter as 
the temperature is lowered is that toward the AFI phase) and below we will explain how 
this can happen, even though U 2 U 3 > 0. 



D. Summary 



Note that we used the amplitude p of the incommensurate wave vector to fix U 2 , in 
contrast to the work of Y, who somehow uses this data to fix Z7 4 . As we have said, the effect 

17 



of U4 is to scale the amplitude of variation of the free energy as one traverses what, when 
U4 is zero, would be the degeneration line. In other words, U4 scales AT, the difference in 
the critical temperatures for FI and AFI fluctuations and the negative sign of U4 is crucial. 
We find that AT is extremely small because it is scaled by the long range interaction U4 
between second neighboring TLL's. 



III. COMPETITION BETWEEN FI AND AFI STATES 

We now analyze the competition between FI ordering (at q z = 0) and AFI ordering (at 
q z = 3tc/c). Although the mean field value of the transition temperature depends only very 
weakly on q z , it is simplest to invoke a model in which only FI fluctuations at q z = and 
AFI fluctuations at q z = 1/2 compete. Therefore we are led to consider the model free 



energy of the form 29] 



To = -(T - T + A/2)x 2 A + ^(T - T - A/2)x 2 F + u[x\ + x 2 F f + vx\x\ , (28) 

where xa {%f) is the AFI (FI) order parameter. The mean-field temperature for AFI (FI) 
ordering is T — A/2 (T + A/2) and A is positive for U4 < 0. Then if, as is usually the 
case, A is temperature independent, one would predict that as the temperature is reduced, 
one would first enter the FI phase, which is not what we want. So we propose a mechanism 
such that A is temperature-dependent, so that as the temperature is decreased, we follow 



the dashed trajectory on the phase diagram for this model 29j shown in Fig. [121 There is 



no reason to expect that within the models considered thus far that A (which arose from 



the value of U4) should have a relatively strong temperature dependence. [30f| It is known 31] 
that the terms of order x 4 in Eq. (128]) implement the fixed length constraint on the variables 
and lead to a temperature dependent renormalization of the coefficients of the quadratic 
terms. But there is no reason to think that such a renormalization will affect the AFI 



order parameter much more than the FI order parameter. It has been suggested 32j that 



this anomalous crossover from FI to AFI fluctuations could be explained by "order from 



disorder." [33(] Here this mechanism (of Ref. |22j) relies on orbital fluctuations. However, it 



would seem that the spin-orbit interaction would cause the orbital degrees of freedom to be 
locked to whatever ordering occurs in the spin degrees of freedom. So we do not consider 
this mechanism. While it is true that at zero temperature quantum fluctuations exist in 
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an antiferromagnet[34J but are zero for a ferromagnet, we are too far from that regime to 



invoke quantum fluctuations. Similar effects do arise from thermal fluctuations. [33| But 
here the antiferromagnetic spin-wave energy is linear whereas the ferromagnetic spin-wave 
energy is quadratic in wave vector. Therefore for identical coupling constants, ferromagnetic 
fluctuations have lower energy than their antiferromagnetic counterparts. This argument 
suggests that FI fluctuations should be stronger than AFI fluctuations. Therefore we reject 



the suggestion 



32| that the cross over from FI to AFI fluctuations can be attributed to this 



mechanism, known as "order from disorder." 
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FIG. 12: (Color online) Mean-field phase diagram 
v > 0. Since the dashed trajectory in the disordered (P) phase is closer to the FI phase than to 
the AFI phase, FI fluctuations dominate AFI fluctuations in the P phase. 



Instead, to obtain the proposed trajectory shown in Fig. H2lwe invoke the coupling of the 
FI and AFI variables to a noncritical variable, Y, so that the free energy is now F = Fq + V , 



where 



35] 



V = a\x z \Y + (1/2)xy 1 Y 2 



(29) 



where here x denotes either the FI or AFI order parameter. Also a is a temperature indepen- 
dent coupling constant and xy is the stiffness associated with Y and is almost temperature 
dependent because Y is far from criticality. (As we shall see, a suitable choice for Y is a 
zone-center phonon.) Since Y is a noncritical variable we can eliminate it by minimizing F 
with respect to it, in which case we obtain 

F = F -(l/2) XY a 2 \x 2 \ 2 . (30) 
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To leading order in the fluctuations we replace \x 2 \ 2 = x(q)* 2 x(q) 2 by 
4x(q)x(q)*(x(q)x(q)*), where (Z) is the thermal expectation value of Z. 31] Then, if 
the coupling constant for FI fluctuations is api and that for AFI fluctuations is ciafi, this 
mechanism leads to the result 

A(T) = A - [a 2 AFI - a| / ] Xy (x(q)x(q)*) , (31) 

where we assume the thermal average is the same for FI and AFI fluctuations. So, if 
dpi <C a afi, then this mechanism leads to a renormalization of the quadratic term which is 
stronger for the AFI fluctuations than for the FI fluctuations. Then the natural temperature 
dependence of the thermal average of |<r 2 | can give the trajectory we desire. 
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Antiferro: f f 

FIG. 13: (Color online) a) The in-plane components of forces on sites L, M, and N in TLLo (as in 
Fig. 1) in the F configuration. The solid line vector forces are from the charges (solid circles) in 
the TLL above TLLo an d the dashed line vector forces are from the charges (dashed circles) in the 
TLL below TLLo- The charges at L, M, and N are negative, positive, and positive, respectively. 
(Red = positive, blue = negative.) The larger dots and heavier lined circles are charges of twice 
the magnitude of the smaller dots and lighter lined circles, b) The net force on the sites assuming 
the separations between all planes are the same. For the AF configuration the dashed forces are 
reversed and the resulting total forces are listed as "Antiferro:". The forces are nonzero at zero 
wave vector only for AF ordering. 

If we choose Y to be a zone-center phonon, the interaction we consider is written schemat- 
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ically as 



V = W Muj d u\ + ■ ^MnDQiQ,] , (32) 



2 



where Mw^ defines the Debye model, u is a phonon displacement, Ufaj) is the heavily 
screened interaction, and Qi is the charge on site i. The factor in Eq. ( 132]) in square 
brackets is the force on site % due to the charge on site j. 

Since we will need the phonon energies, we implemented a a first-principles calculation of 
the energies of the zone-center phonons in LuE^O^ The calculations were performed within 
the plane-wave implementation of the local density approximation to density functional the- 
ory as implemented in the PWscf package. [36j We used Vanderbilt-type ultrasoft potentials 
with Perdew-Zunger exchange correlation. A cutoff energy of 408 eV and a 9x9x9 k-point 
mesh were found to be enough for the total energy to converge within 0.5meV/atom. The 
gamma phonon energies were calculated using the supercell method with finite difference. 37] 
The primitive cell was used and the full dynamical matrix was obtained from a total of eight 
symmetry-independent atomic displacements (±0.02 A). 

The primitive cell contains one formula unit of LuFe204, giving rise to a total of 21 phonon 
branches. The phonon modes at T are classified as T (q=0) = 4A 2u (IR) + 3Ai g (K) + 4E U (IR) 
+ 3E g (R), where R and IR correspond to Raman- and infrared-active, respectively. The 
nondegenerate (A) and doubly degenerate (E) modes correspond to motion along the c-axis 
and within the ab-plane, respectively. In the Raman-active modes the atoms at (0, 0, z) and 
(0, 0, —z) move out-of-phase (i.e. opposite) while in the IR-modes they move in phase. 

The calculated mode energies and their symmetry labels are listed in Table [III We hope 
that our calculations will initiate more experimental work such as Raman/IR measurements 
to confirm the T phonon energies that we calculated here. From this table we see that the 
characteristic phonon energies are of the order of 50 or more meV. Since the coupling to the 
phonons with the lowest energies will be the most effective, we show the two lowest energy 
phonons schematically in Fig. [T41 

The lowest energy mode with E u symmetry and 11.4 meV energy corresponds to displace- 
ments in the ab plane in which the LuO and FeO-bilayer move in opposite directions as rigid 
units (see Fig. [Ilk). Hence the energy of this mode is basically determined by the strength 
of LuO-Fe-0 bond angle. Even though this mode has the lowest energy, its symmetry is 
not right to create the electrostatic forces needed for our mechanism. The next mode has 
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(a) E u 11.41 meV (b) E 19.99 meV 




FIG. 14: (Color online) The displacements of the two lowest energy modes as discussed in the 
text. Here O ions are red, Fe ions are purple, and Lu ions are green. The light blue ovals indicate 
units which move approximately rigidly. 

the Eg symmetry and it corresponds to displacements in the ab-plane (Fig. H4"b) in which 
the two TTL's of each bilayer moves in opposite directions while LuO-layer is fixed. Hence 
this modes involves twice as much O-Fe-0 bond bending as the lowest energy mode and 
interestingly it has about the twice energy (20 meV) of the E u mode shown in Fig. [l~4"a . As 
we shall see, it is this mode that creates the electrostatic forces needed for our mechanism. 

TABLE II: List of phonon symmetries and calculated energies (in meV) of LuFe204 at the T point 
of the primitive cell, as obtained from the first-principles calculations described in the text. 



Mode Symmety 
Energy (meV) 


E u Eg Ai 9 E u 
11.41 19.99 20.02 31.48 38.46 41.20 


Mode Symmety 
Energy (meV) 


Alg Eg E U Eg A\g 

53.25 54.73 57.69 58.79 62.62 84.96 



Accordingly, we look for charge-phonon interactions which involve zone-center transverse 
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(to the c-axis) phonons. We now analyze the force on the Fe charges in a given TLL, which 
we denote TLLo, from the nearest neighboring Fe TLL's above and below TLLo. Since we 
consider coupling to the lowest energy phonons, which involve motion transverse to the c 
axis, we will only consider forces in the plane of the TLL. Phonon modes which decrease 
the distance between TLL's will involve higher energy. One sees that a low energy mode 
which can couple the way we want is a zone-center phonon in which alternate TLL's are 
displaced transversely relative to one another. For this rhombohedral lattice, such a mode 
is an E g mode. As mentioned this is the E g mode at 19.99 meV. To see whether this mode 
couples differently to FI and AFI fluctuations, we have only to analyze the transverse force 
on one TLL from the TLL's above and below it. Since we are near the transition, we assume 
the R3 structure of the fluctuations (see Fig. El where we choose = 0) and add up the 
forces in Fig. [131 Also, we simplify the argument by neglecting the fact that the interlayer 
separation is different for TLL's within the same bilayer and for TLL's in adjacent bilayers. 

We now estimate quantitatively the effect of this coupling in Eq. ( 1321) . As in Fig. 1, 
the Q's are given in terms of the order parameter x%, where X indicates either the FI or 
AFI configuration of TLL's. Because the transverse motion of planes is relatively soft, we 
consider displacements to lie within the TLL plane. When minimized with respect to Ui, 
the free energy is 

4M 4 (U\ 2 fr^ 2 



= K7) Ki) « • (33) 

where £ = (r/U)(dU/dr), rn is the component of r within the TLL, and z is the effective 
number of nearest neighbors. Also we set U = (U2 + U%)/2 = 16. 5K pa 1.5 meV. The actual 
number of nearest neighbors is 6, but since the forces do not all add up, we take z = 3 for the 
AF configuration and z = for the F configuration where the forces from adjacent TLL's 
nearly cancel. We set Hojd = 20 meV, r = 5A, and rn = 1°A. The factor £ depends on how 
rapidly the interaction decreases with distance. For bare Coulomb interactions £ = — 1. But 
we are far from that regime. We take £ = —10, which is a value often found for exchange 



381 ] Also the Fe mass is M = 60 amu, so its reduced mass is 30 



interactions in insulators, 
amu. Thus 

F Ph = -F |4 F | ->• -^(8F (\x AF \ 2 ))\x AF \ 2 , (34) 
where F pa 0.004meV. Then identifying {8F {\xaf\ 2 )) as the renormalization of Tqq we get 
A(T) w A - 0.032meV(|4 F |) w 0.04K - (0.4K)(|4 F |) . (35) 
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which is enough to shift ordering from F to AF at Tqo, where (I^afI) ~ 1/2. 

Finally, we should mention that this frozen phonon occurs whether or not the CO phase 
is commensurate because its origin is in a coupling of the form 

V ~«(q = 0)<7(q)<r(-q) , (36) 

where cr(q) is the CO order parameter. 



IV. THE MAGNETIC PHASE TRANSITION 



A. Phase of the R3 Structure 



We now discuss the magnetic phase transition at which SO appears. At first we neglect 
the fact that the system is a mixture of spins of magnitude 2 and spins of magnitude 5/2 
and we assume that the uniaxial anisotropy aligns the spins along the c axis. Then one 
introduces the local order parameter S(r) as the thermal average of S z (r), the ^-component 



of spin at the site r. Also, as noted in Ref. [ll 



, if one neglects the coupling between spin and 



charge, the symmetry of the SO free energy is the same that of the CO free energy of Eq. 
(2). Thus, if the transition is assumed to be continuous, the ordering wave vector for this 
transition should be unstable relative to the X point, just as we argued (in connection with 
Fig . 6) in the case of the CO transition. In that case the representation analysis of Ref. 



11 



for the wave vector (1/3, 1/3,0) is not definitive. However, we temporarily overlook the 
possible instability of the X wave vector and apply Landau theory to the phase transition 
as if this wave vector were stable. (In Appendix A we discuss some difficulties in applying 
representation analysis to this transition.) Therefore we write the free energy in terms of 
Si(q) = 5*2 (q) = S(q), where the subscript labels the two Fourier components of the unit 
cell. We have that 

F = (1/2)(T - T S o)\S(q)\ 2 + u\S( q )\ 4 + t;|S(q)| 6 + . . . + wS^f + ^(q)* 6 , (37) 

where Tso is the magnetic (SO) transition temperature and we will set 

5(q) = |S(q)|e*. (38) 

Under inversion symmetry S(q) — > 5*(q)*, which implies that w in Eq. (1371) is real. The last 
term in Eq. fl37|) is the lowest order one that fixes the phase <fi of the order parameter. (It 

24 



should be noted that it is not easy to fix the this phase using only scattering data.) There 



are two cases: 



(p = nn/3 ii w < 



(n + l/2)7r/3 if w > ; 



(39) 



with the results for the amplitudes in the magnetic unit cell as given in the caption to Fig. 
5. To determine the net moment of these structures it is necessary to analyze the admixture 
of wave vector (0, 0, 0).|40j Such an admixture comes about because q = (1/3, 1/3, 0), is 1/3 
of a reciprocal lattice vector and this fact allows an additional term, V, in the free energy, 
where 

V = -5(0, 0, 0) [aS(l/3, 1/3, 0) 3 + a*S(l/3, 1/3, 0)* 3 ] + ^^5(0, 0, 0) 2 , (40) 

where x is a stiffness (which is nearly temperature independent) and a is a constant which 
must be real in view of inversion symmetry. Minimizing V with respect to 5(0, 0, 0) we find 
that 



5(0, 0, 0) = a X [S(l/3, 1/3, 0) 3 + 5(1/3, 1/3, 0)* 



(41) 



If w > 0, then Eq. ( l39i) indicates that 5(0, 0, 0) is zero, whereas for w < the system has 
a nonzero net moment, M. The early data of Ref. Ig suggests that M is nonzero. However, 
recently we have learned 41] that the system is more likely to have M = 0, in which case we 
must choose w > 0. In this case one of the three sublattices is disordered. (See the caption 
to Fig. [5] with cf) = tt/6.) In this structure, all spins within a plane perpendicular to the 
ordering wave vector have the same value, a, —a, or 0. This type of partial ordering was 
observed in the orientational ordering of solid methane 42j, and, as in that case, we would 
not expect a phase with partial disorder to continue to exist to arbitrarily low temperature. 
In Appendix B we obtain M(H) for this structure. 

If, instead, the case w < is realized, then one would have 



5(0,0,0) = 2a X S(l/3, 1/3, Of 



and if 



5(1/3, 1/3,0) = B~(T so -T)P 
then, within mean field theory 

5(0,0,0)^A~(Tso-T) 3 ^, 



(42) 



(43) 



(44) 
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which gives an unusually large critical exponent for the magnetization. For liquid 



crystals, 



43| this effect has been analyzed in detail within the renormalization group. An 



effect similar to this has been seen for CO. 
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Next we discuss the magnetic eigenvector which was chosen in Ref. [ll to be [1, 1] (to 
best fit the experimental data). With this choice of eigenvector the spins form planes 
(perpendicular to the wave vector) of spins with amplitudes proportional to 0, —1, and 1. 
How is this choice to be justified within Landau theory? In the 'minimal' model used for 
CO, one sees that for p = (the commensurate case) Fu=0 and one has isotropy in S\ z , 
5*2,2 space so that the eigenvector can be [cos 6, sin 8] with any choice for 9. One way to 
explain that the eigenvector is [1, 1] is to invoke an interaction which tends to make the two 
spins in the rhombohedral unit cell parallel, so that _F 12 is negative real. Since the spins are 
aligned along z, the dipole interaction could accomplish this. However, the energy V of this 
interaction is probably too small: 

V = = _ .045K , (45) 

where we set g — 2, S — 5/2, and r = 2d' + d = 14. 1A and combine with a second term 
for which r = d! + 2d = 11.1 A. Alternative mechanisms to stabilize the antiferromagnetic 
spin structure involve thermal fluctuations or (as we discuss below) the distortion due to 
the frozen E g phonon. 

Finally we discuss the diffraction at half integer values of L which has a magnetic signature 
and which requires positing a magnetic unit cell which is doubled along the c direction. As 



stated in Ref. 



11 



one can think of this additional diffraction as being due to "the charge 
ordering, which decorates the lattice with differing magnetic moment on the Fe 2+ and Fe 3+ 
sites..." This effect can be seen within Landau theory as follows. We introduce an additional 
free energy V of the following form, consistent with symmetry, 

V = -a]>>(r)^(r) 2 , (46) 

r 

where a is a constant. The effect of this term is to increase (decrease) ^(r) 2 when the site 
is occupied by an Fe 3+ (Fe 2+ ) ion. In Fourier transform language this is 

V = -a £ x n (q)5 n (k 1 )5„(k 2 )A(q,k 1 ,k 2 ) , (47) 

n;q,ki,k2 

where A enforces wave vector conservation modulo a reciprocal lattice vector and the sub- 
scripts label the Fe sublattices. The term we focus on here involves charge ordering, so that 
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q = (1/3, 1/3, 1/2) + 8, where 8 is the incommensurability. The magnetic variables then 
can involve the wave vectors ki = (1/3,1/3,-1/2) — 8 and k 2 = (1/3,1/3,0). Then we 
see that this interaction couples S(l/3, 1/3, 0) and 5(1/3, 1/3, —1/2)* + 8. Thus the critical 
magnetic eigenvector is a mixture of these two variables. This corresponds exactly to the 
idea of "decoration," but it is hard to estimate the importance of this effect. 



B. Removal of Frustration 



To develop further insight into this frustrated spin system, it is useful to recall the results 
for the magnetic structure of the rhombohedral /3-phase of solid oxygen whose structure only 
differs from LFO in that there are no bilayers: all TLL's are equally spaced. (For a review see 



Rcf. 
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A convincing theoretical analysis based on quantum spin-wave theory was given in 



Refs. |18| and[l9|. However, when various theoretical results were experimentally tested [4 71]. 

it was not entirely clear which theoretical model was most appropriate for /3-oxygen. In any 

I I 

event, the magnetic correlation length is so short (5A)[47J that is seems unrealistic to speak 
of any long range order. 

Accordingly, a central open question is to explain why the SO in LFO is so different from 
that of /3-oxygen. One possibility is that the small distortion which we invoked to explain 
the crossover from FI to AFI might introduce small addition exchange interactions which 
remove the frustration of the rhombohedral antiferromagnet. To explore this possibility we 
write 

J n = J<® + 6 J n , (48) 

where the J^, s are the exchange interactions which are consistent with the R3m symmetry. 
In Fig. [15] we show a set of interactions 5J n which have the correct symmetry to be induced 
by the frozen E g phonon and which, if they are dominant, resolve the frustration. Note that 
even though CO is incommensurate, the frozen phonon is commensurate. So here we are 
considering a commensurate effect of incommensurate charge ordering. 

To analyze this possibility, we replace the J3 interaction by Jz(l + e) for the dashed bonds 
in Fig. [15] and the other J 3 interactions by J^(l — e). Similarly, we replace the J 2 interaction 
by ^(l — e) for the solid bonds in Fig. [15] and the other J 2 interactions by J 2 (l + e). Then 
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a b c 

z 

FIG. 15: (Color online) Additional exchange interactions attributed to the transverse frozen 
phonon. Dashed (full) lines denote additional antiferromagnetic (ferromagnetic) interactions. The 
values of S z for all sites in the same y-z plane are a, b, and c as indicated, where b = —a and c = 0, 
in the notation of Eqs. (I42p and (I43D . It is assumed that the in-plane interactions are dominant. 

Eq. fllj) remains valid but now 

A = (1 + e y°4y^/e co ^ aqx /2) + (1 - e ) e -^v / 3/3 

= -2e - (V3/2)( Px + ipy) - e Px V3/2 + tep y V3/6 + 0(p 2 ) . (49) 

To maxmimize |F 2 i| (for U2U3 > 0) set exp(icq z /3) = A/A*, so that 

\F 21 \ = (U 2 + U 3 )\A\ . (50) 

Then the minimal eigenvalue is found by minimizing 

p(p x , Py ) = T-3f/ 1 + (3/4)f/ 1 p 2 

~(U* + U 3 ) ([2 e + (1 + 6)(v / 3/2)p,] 2 + (3/4)(l - e /3) 2 p^ . (51) 

When e = 0, this is a function of p x +Py and is consistent with the existence of a degeneration 

28 



line. However, when e is nonzero, then the minimum occurs for p y 

p x = (l + e)(U 2 + U 3 )/(V3U 1 ) . 



= (so that q z 



= 0) and 
(52) 



One might wonder where the extrema we found for CO at cq z = 2pir } for p = 1 and p = 2 
have gone. The answer is that there are three CO domains corresponding to which there are 
three distortions, the transverse phonon displacement being perpendicular to the in-plane 
projection of the ordering wave vector. So we have three different SO domains, each one 
tied to one of the three possible CO domains. 

An important question is: since the wave vector is not at a special, high-symmetry 
point (see Fig. 6), the wave vector should not be commensurate if the SO transition is a 
continuous one. It is not clear that the sensitivity of the neutron diffraction experiment 
of Ref. llll is sufficient to detect the very small incommensurability that might attend this 
magnetic transition. It would be of interest to have a high precision determination of the 
SO wave vector, to check whether it is or is not commensurate. 

If the SO phase is truly commensurate, then one would have to entertain a scenario to 
accommodate such a fact. The one scenario that is excluded is that the commensurate 
SO state is reached via a single continuous phase transition. Possibly there are two phase 
transitions, the first, in which there develops incommensurate order, followed by a second one 
into the commensurate antiferromagnetic state. 45J The presence of two nearby transitions 
in parameter space would seem to signal a nearby multicritical point at which the two 
transitions coincide. To check for that, it would be useful to have very precise measurements 
of the specific heat and susceptibility to get the critical exponents that characterize this 
transition. A different scenario is that the magnetic transition is a first order one to a 
commensurate SO phase. 



C. Field Cooling 
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where it is shown that 



Finally, we mention the intriguing data of Fig. 3c of Ref. 
cooling in a field along (1, —1,0) from T > Tco causes a pronounced reduction in the CO 
diffraction at 300K. This data raises a natural unanswered question: does this 'missing" 
intensity in the AFI scattering show up as new FI scattering at (1/3, 1/3, n) for integer n. 
If so, it would mean that the magnetic field could stabilize the FI state for T < Tco an d it 
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would be of interest to know whether such a state was or was not commensurate. 

Here we present apartial explanation for the above field cooling scenario. We start by 
noting that in Ref. [l3j it is shown that application of a magnetic field in the plane of the 
TLL's decreases the AFI correlation length. This suggests that such a magnetic field tends 
to destabilize the AFI phase, possibly making the FI phase relatively more stable. If this 
is the case, the one might have a phase diagram like that shown in Fig. [161 Then in the 
various scenarios of cooling one would start in the disordered phase at points like A and B 
and cool to points A' and B'. Clearly, if this is done at zero field and then a field not large 
enough to go into the FI phase is applied, no dramatic field dependence will be observed, in 
agreement with their observations. However, if one starts from point like A or B, then when 
the system is cooled it passes though a region of FI ordering which then can be supercooled 
while reaching the final points A' or B'. Then, as a function of time the system would 
evolve in some irregular process into the equilibrium state of AFI order. This might happen 
without displaying a dramatic dependence on the field-cooled value of the magnetic field. 
Indeed the data shows that after a sharp decrease for very small field, the resulting AFI 
order does not depend strongly on H. This proposal suggests that it would be useful to 
monitor the (1/3,1/3,0) reflection under field cooled conditions (to see if the decrease in 
intensity at (1/3, 1/3,3/2) is accompanied by an increase at (1/3, 1/3,0). It would indeed 
be interesting if an in-plane magnetic field could stabilize a nonzero polarization. Note that 
an in-plane magnetic field may be more effective than one parallel to the c-axis, because the 
perpendicular susceptibility is usually larger than the parallel susceptibility. 

V. SUMMARY 

Here we briefly summarize our conclusions. 

• We show that the appearance of an incommensurate wave vector for charge ordering is 
a result of symmetry (or more accurately, due to a lack of symmetry). 

• By comparing our theory with experiment we have assigned values to several of the 
phenomenological charge-charge interactions. In particular, the signs and magnitude of the 
incommensuration and the fact that ferroelectric fluctuations dominate in the disordered 
phase are explained by simple choice of these interactions. 

• The cross over, as the temperature is lowered through the charge ordering temperature, 
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FIG. 16: (Color online) Proposed phase diagram for CO as a function of maganetic field H and 
temperature T. Dashed lines represent continuous transitions and the solid line a first order 
transition. 



from ferroelectric to antiferroelectric incommensurate structure can be explained by the 
temperature dependent renormalization of the transition temperature due to charge-phonon 
coupling. 

• We have performed a first-principles calculation of the zone center phonon energies 
(assuming no ordering of charge or spin). The phonon with the correct symmetry to couple 
effectively to the charge ordering has a rather low (20 meV) energy, corresponding to the 
sliding (transverse to the c axis) of alternate Fe layers with respect to one another. 

• We have developed a Landau theory which describes the phase of the recently observed 
spin ordered state having zero net magnetic moment. 

• In principle, if the spin ordering transition is continuous, the spin ordered phase would 
be expected to be incommensurate. This suggests the need for a high precision determination 
of the spin ordering wave vector to check whether it is or is not commensurate. If the spin 
ordered phase truly is commensurate, then it would be of interest to investigate the scenario 
of ordering, which can not be via a single continuous transition. One way to pin down the 
scenario would be to determine the critical indices a, (3, and 7, associated respectively with 
the specific heat, the magnetic order parameter, and the susceptibility. 
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• We have also suggested experiments to test our proposal that the sharp decrease in 
intensity of antiferroelectric charge scattering as a function of magnetic field in the field 
cooled scenario might indicate that the magnetic field tends to stabilize ferroelectric charge 
ordering and possibly a consequent polarization. 
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Appendix A: Representation Theory 



In Refs. 



11 



and 



12 



representation theory is used to analyze possible magnetic ordering 
patterns and charge ordering patterns, respectively. In their approach, they implicitly as- 
sume that the wave vector at the appropriate X point is stable with respect to the addition 
of further neighbor interactions. As we have seen in Sec. II, this assumption is not actually 
valid, especially for CO. To see this explicitly, consider the structure of the two by two 
matrix F nm of Eq. (4) which determines the eigenvectors. Exactly at the X point and when 
arbitrary interactions are allowed, F 2 \ is scaled by the interactions between sites ^1 and 
j^l which are displaced from one another by a vector along the c-axis. This interaction is 
extremely small, since it connects sites which are not in adjacent bilayers, but are in second 
(or further) neighboring bilayers. Representation theory bases the eigenvector equation on 
this symmetry and leads to eigenvectors that are either even or odd under inversion. 

However, as noted in Sec. II, this type of analysis is invalidated by the fact that for 
LFO the X point is not actually stable. For wave vectors near the X point, we explicitly 
displayed in Eq. (7) the term in F 12 which is linear in the displacement from the X point. 
The interactions U2 and U3 which scale this linear term are very much larger than V\ t 2 
whose existence is ignored if representation theory is invoked for the commensurate case. 
The major effect of this linear term is that the eigenvectors, instead of being even and odd, 



as in the analysis of Refs. 



11 



and 



12 



are now complex and are determined by the phase of 
F21 given by Eqs. (4) and (7). Note that this phase will, in general, be different for each 
of the three domains, and inclusion of the correct phases, might affect the determination of 
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the domain populations. 



Appendix B: Equation of State for the Antiferromagnet 

In this appendix we obtain M(H), where M is the net magnetization (along the c axis) for 
the model of Eq. (|37|) and H is the external field applied parallel to the c axis. Accordingly 
we add to the free energy the term —HM, where x is the parallel susceptibility, so that with 
|if>(q)| = cr, we have 



F= -{T - T so )a 2 + «<t 4 + Awa 6 cos 2 (30) - 2aa 3 M cos(30) + -x _1 M 2 - ME , (Bl) 



1 

2< 



where x is the parallel susceptiblity and we kept ^-independent terms only up to order a 4 
because our analysis is not valid when T <C T$o- Minimizing with respect to M yields 



M = x[H + 2aa 3 cos(30)] 



(B2) 



so that 

F= -(T -T so )a 2 + ua* + 4wa 6 cos 2 {3<f)) --x[H + 2aa 3 cos{3<f))} 2 . (B3) 
2 2 

When this is minimized with respect to <fi we find two regimes: 



H>H a 



H<H C : 



M = x 



H + 2aa c 



where 



M = X H 1 + 



{Aw - 2a 2 x)cr 3 
ax 



2a 2 X 



Aw — 2a?x 



(B4) 



(B5) 



which leads to the M versus H curve shown in Fig. [T71 The value of a is approximately 
[(T so - T)/(4 M )] 1 /2 f or T near T so . 
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FIG. 17: (Color online) The c-component of magnetization M (solid line) versus the magnetic field 
H along c for the antiferromagnetic phase. Here Mq = 2ax& 3 and H c is given by Eq. (|B5p . Near 
T so both M and H c are of order (T so - T) 3 / 2 . 
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